Stochastic number projection method in the pairing-force problem 
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A new stochastic number projection method is proposed. 
The component of the BCS wave function corresponding to 
the right number of particles is obtained by means of a 
Metropolis algorithm in which the weight functions are con- 
structed from the single-particle occupation probability. Ei- 
ther standard BCS or Lipkin-Nogami probability distribu- 
tions can be used, thus the method is applicable for any pair- 
ing strength. The accuracy of the method is tested in the 
computation of pairing energies of model and real systems. 

PACS number(s): 21.60.-n, 21.10.Ma, 02.90.+P 



The microscopic model of superconductivity intrcu 
duced in 1957 by Bardeen, Cooper, and Schrieffer (BCS)El 
has had astonishing success in correlating and explaining 
the properties of simple superconductors in terms of a few 
experimental parameters. For its conceptual and compu- 
tational simplicity, the BCS method has been widely used 
as the first step ia nuclear structure calculations involv- 
ing pairing forcesa. The theory is quite satisfactory when 
the number of valence nucleons is large and the pairing 
interaction is strong (as compared with the level spac- 
ing). But in a nucleus with a relatively small number of 
valence nucleons, or with a large spacing between levels, 
the BCS method fails. As it is well known, the method 
has two inherent drawbacks: 

1) The BCS wave function is not an eigenstate of the 
number operator. The energy obtained from this wave 
function is , therefore, biased with an inaccuracy caused 
by the number fluctuation. 

2) In some cases, there may be a critical value of the 
pairing force below which the BCS equations have no 
non-trivial solution. JUxact calculations show that this 
behaviour is spuriousu. 

Attempts have been made at improving the method. 
Kerman et alu showed that the wave function obtained 
by projecting out the BCS function to the sector with the 
right number of particles is a very good approximation to 
the exact wave function. Many works have been devoted 
to particle-number conserving approaches, including pro- 
jection before and after variationcTtll, and many others, 
see for example Ref. [l2| However, it is well known that 
these approximate methods can lead, in some cases, to 
significant errors. 

The exact solution of this problem is only available for 
some very simple systems, such as a single level or a two- 
level modeSEj, or the caaenof equispaced doubly degep. 
erate single-particle levelst 3 ! Richardson and Shermantil 



have also developed a general method for determining 
the exact eigenvalues and eigenstates of the pairing force 
Hamiltonian when the pairing strength is constant. In 
a recent paper, Cerfta has proposed the application of a 
diffusion Monte Carlo (DMC) technique to the treatment 
of the pairing force in nuclei. By making use of DMC, it 
is possible to compute exactly the ground-state energy of 
the system with a general, state-dependent, pairing inter- 
action, at least in principle. However, the DMC method 
can only deal with pairing Hamiltonians in which the in- 
teraction matrix elements are strictly positive. This is 
always the case for the nuclear pairing Hamiltonian, but 
other types of pairing problems, involving Coulomb in- 
teracting particles for example, can not be solved in this 
way. 

In this paper, we present a novel-stochastic algorithm, 
based on the Metropolis methodlia, for projecting out 
the component of the BCS wave function with the cor- 
rect number of particles. The method does not depend 
on the type of interaction involved, as long as we can 
assume that the particular BCS functional form holds 
for our ground-state wave function. Moreover, number 
projection canJjifi done starting either from the BCS or 
Lipkin-NogamiQij calculated occupation probability dis- 
tribution. Therefore, it can be applied for any strength 
of the pairing interaction. In what follows, we first in- 
troduce the model Hamiltonian, then discuss our Monte 
Carlo projection method. Next, we apply the method 
to simple systems, and compare the results with those 
obtained from other approaches. 

We consider the problem of a many-body system de- 
scribed by the Hamiltonian, 
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The at and at create particles in time-reversed con- 
jugate single-particle states \j) and \j) with energies tj. 
The interaction V pa i r scatters only time-reversed pairs of 
particles from the occupied levels \ j, j) to the empty ones 
\j' , j'). The indices j and j' run from 1 to fi, where Q, is 
the total number of conjugate orbit pairs. Our purpose 
is to compute the exact N-particle ground-state energy 
of H. We construct the ground state of H by projecting 
BCS-type wave-functions onto the N-particle sector. Let 
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us define the number of pairs, N p = N/2. The following 
projected function is obtained: 
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are the normal- 



isation coefficients, |Vac) is the vacuum state, Vj is the 
amplitude for finding a pair of particles in time-reversed 
levels, and Uj-the amplitude for the levels being 

empty. These amplitudes are obtained either by BCSEl or 
Lipkin-Nogami (LN)ltU methods. Using the Hamiltonian 
given by Eq. ([!]) and the wave function (|^), we arrive to 
the following expression for the projected ground-state 
energy, 
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W(ji,...,j Np ) s{ji,...,j Np ), (3) 



where the sum runs over all possible combinations of 
N p occupied single-particle double degenerated states 
{ji , . . . , j jv p } , from a maximum of O states allowed in the 
Monte Carlo evaluation. The "weights", W, and "ener- 
gies", e(jx, ...,Jn p ), are defined as 
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This summation is impracticable as it involves a huge 
number of terms, ( jv ) > f° r usual model spaces. For- 
tunately, the expression (|^) for the projected energy 
allows a simple Monte Carlo evaluation, where the 
ensembles {j% , . . . , jjsr } are generated with probability 
W(ji, . . . , Jn p ) by means of a Metropolis algorithmll3. 
Other equivalent forms of Eq. (||) , see for example Ref. 
[l?], are not suited for the Monte Carlo evaluation. 

In the Metropolis evaluation, we start from the un- 
perturbed ground state and perform 10 5 thermaliza- 
tion steps. We define a trial Metropolis move as a 
random transition of one pair from the occupied lev- 
els \j,j) to the empty ones \j',j')- The trial move is 
accepted or rejected according to the Metropolis rule 
w tii, ■ ■ ■ J'n p )/ w Ui, ■ ■ ■ ,3N P ) > 1, where 7 is a ran- 
dom number uniformed distributed between and 1 . An 
acceptance ratio R is obtained from the thermalization 



loop. Using this ratio a total number of Metropolis steps, 
NMetropoiis = 300 N ' decor relation is used to estimate av- 
erage values of ground-state energies. The number of 
decorrelation steps were taken as: N^eeorr elation — 50/ R. 

In the following, we will compare results for jSeyeral 
soluble, models with those obtained by BCS, LNLjEI and 
DMCE-3 methods. We are interested in the two-body pair- 
ing energy, E pa i r , which we define as E pa i r — E( 
Erv__ ir jtQ), where E( Vpr ,^m and E/ 



'-(v pa , r =o) ~ 

'^(v patr ^o): wnere -C(y paii . =0 ) ™i ^(y pair ^ ) ar e, respec- 
tively, the ground-state energies of the system without 
and with pairing interaction. 



A. The symmetric two-level model 

In this paragraph, we study an exactly solvable, sym- 
metric, two-level model, with number of particles, TV, 
and level degeneracy, fl. The pairing interaction is taken 
constant and equal-to G. This model has been first exam- 
ined by Hogaasenlla, and its exact, solution was studied 
in detail by Rho and Rassmussera in the case N = Q. 
More recently, the general case N ^ f2 has been dis- 
cussed in Rcf. [h]. The exact pairing energy E pa i r is ob- 
tained by introducing two sets of quasispin operators, 
so that the problem finally preduces to the diagonaliza- 
tion of a tri-diagonal matrixLniy. The results of the BCS 
approximation, the LN prescription and the DMC and 
Metropolis projection methods (MCP) using Vj and Uj 
either from BCS (MCP BCS ) or LN {MCP LNM ) calcula- 
tions are compared against exact results in Table [| for 
O = 10 and N = 4, 10 for several values of the interac- 
tion strength, G. It is shown that Metropolis Projection 
methods give very good agreement with the exact results 
(the same quality of agreement as the DMC method). 
Recall that for N t^-i^ the BCS ansatz always have a 
nontrivial minimumtS. For N — f2, a nontrivial solution 
is found only for G greater than 1/(N — 1). 



B. Equidistant doubly degenerate levels 

In this case, we deal with a system of equispaced 
doubly degenerate single-particle levels and a constant 
pairingj—ipteraction. This problem has been solved 
exactlyLJO for some model spaces with values of the 
interaction strength reproducing typical nuclear pairing 
energies. Results obtained from different methods are 
shown in Table || for SI = N = 8 and three values of G. 
Again MCP is in good agreement with the exact results. 



C. The 100 Zr nucleus 

In this section, we apply our method to compute the 
ground-state pairing energy of the strongly deformed, 
neutron-rich nucleus 100 Zr. Recently, the pairing prob- 
lem in this nucleus was studied by the DMC methodEI 
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The average field was assumed, to be an axially de- 
formed Woods-Saxon .pptentialEj with Cassinian ovals 
shape parametrizatior£3. The universal Woods-Saxon 
parameters proposed by Dudek et alEij were used in the 
single-particle level calculations, except that smaller val- 
ues (i?o = 1-25 for both particles) ofr-the central poten- 
tial radius par.ameter were employedcJ. The computer 
code CASSINo using Cassinian ovals parametrization 
was applied to obtain single-particle level energies. In 
the pairing calculations, we keep all single-particle levels 
(f2 = 20) obtained in the deformed ground-state poten- 
tial in a 10 MeV interval around the Fermi level (5 MeV 
above and 5 MeV below). The single-particle energies 
relative to Fermi level energy are listed in Table [II . The 
number of pairs in this system is equal to the number of 
considered single-particle levels below the Fermi energy, 
i.e. 9 (number of particles N — 18). The BCS , LN and 
MCP results are c omp ared against benchmark DMC cal- 
culations in Table IV for interaction strength G — 0.255. 
Differences between publishedE2l and present BCS results 
arise from the consideration of the self-energy. The ac- 
curacy of MCP is, once more, comparable to DMC. 

In conclusion, we have proposed a method for project- 
ing out BCS-like wave functions to the A^-particle Hilbert 
space. The method shares with the BCS approach its 
computational simplicity. The use of LN approximate 
projection makes it possible to apply the MCP at any 
pairing strength, even below the BCS critical coupling, 
if it exists. These properties make the method very suit- 
able for nuclear structure calculations. We stress that, 
unlike DMC, any pairing interaction can be treated in 
our approach. 

Besides the nuclear pairing problem, there are other 
possibilities of application of the MCP method . Re- 
cently, we have used it in the computation of the ground- 
state energy of electron-hole systems in a quantum dotc3. 
The electron-hole attractive Coulomb potential is the 
pairing interaction in this case. MCP is used to im- 
prove the BCS estimation. Th e BCS approach to electron 
correlations in moleculescllca is another example with 
Coulomb matrix elements. It was shown that this ap- 
proach fails to reproduce the correlation energy of small 
molecules, probably because of the need of an exact pro- 
jection. MCP may be a good alternative. Finally, we 
shall mention as a possibility of application the recent 
study of superconductivity in ultrasmall grainsEH 

Some of these problems are currently under investiga- 
tion. 
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TABLE I. Pairing energies for selected values of the inter- 
action strength in the two-level model with Q — 10 double 
degenerate levels and N = 4, 10 particles. 



G 


Exact 


BCS 


LNM 


DMC 


MCPscs 


MCPijvM 


n = 10, N = 4 


0.024 


0.202 


0.161 


0.205 


0.205(3) 


0.203(3) 


0.205(3) 


0.138 


1.495 


1.269 


1.507 


1.498(3) 


1.500(5) 


1.500(5) 


1.105 


18.05 


16.276 


18.05 


18.05(9) 


18.05(2) 


18.05(2) 


Q = 10, iV = 10 


0.024 


0.147 


0.000 


0.122 


0.146(2) 




0.145(5) 


0.060 


0.36 


0.000 


0.311 


0.373(9) 




0.353(9) 


0.100 


0.711 


0.000 


0.603 


0.717(5) 




0.696(7) 


0.500 


10.56 


9.300 


10.56 


10.56(5) 


10.56(5) 


10.56(5) 


1.105 


28.40 


25.64 


28.40 


28.40(7) 


28.40(2) 


28.40(2) 



TABLE II. Pairing energies for selected values of the inter- 
action strength in the system with Q — 8 double degenerate 
levels and N = 8 particles. 



G 


Exact 


BCS 


LNM 


DMC 


MCPscs 


MCPijvM 


0.7 


5.309 


3.89 


5.245 


5.27(3) 


5.28(5) 


5.25(6) 


0.9 


8.018 


6.18 


7.977 


8.06(2) 


7.98(5) 


7.97(4) 


1.1 


11.181 


8.83 


11.05 


11.05(4) 


11.11(5) 


11.09(9) 



TABLE III. Single-particle energies of 100 Zr in a 10 MeV 
interval around the Fermi energy. The deformation parameter 
is e = 0.33. 



Levels above Fermi energy 


Levels below Fermi energy 


4.684 


0.0000 


4.559 


-0.254 


4.275 


-1.303 


3.566 


-1.520 


3.300 


-2.224 


2.908 


-2.826 


2.235 


-3.914 


1.928 


-4.125 


1.497 


-4.605 


1.323 




0.311 





TABLE IV. Ground-state pairing energy for the interac- 
tion strength G = 0.255 in the 100 Zr nucleus. 



G 


DMC 


BCS 


LNM 


MCPscs 


MCPijvM 


0.255 


5.0 


3.87 


4.94 


5.1(1) 


4.98(6) 
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